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We propose a new approach to probing ergodicity and its breakdown in quantum many-body sys¬ 
tems based on their response to a local perturbation. We study the distribution of matrix elements 
of a local operator between the system’s eigenstates, finding a qualitatively different behaviour in 
the many-body localized (MBL) and ergodic phases. To characterize how strongly a local perturba¬ 
tion modifies the eigenstates, we introduce the parameter Q(L) = (In (V n m/S)}, which represents a 
disorder-averaged ratio of a typical matrix element of a local operator V to the energy level spacing, 
S', this parameter is reminiscent of the Thouless conductance in the single-particle localization. We 
show that the parameter G(L) decreases with system size L in the MBL phase, and grows in the 
ergodic phase. We surmise that the delocalization transition occurs when G(L) is independent of 
system size, G(L) — Q c ^ 1. We illustrate our approach by studying the many-body localization 
transition and resolving the many-body mobility edge in a disordered ID XXZ spin-1/2 chain us¬ 
ing exact diagonalization and time-evolving block decimation methods. Our criterion for the MBL 
transition gives insights into microscopic details of transition. Its direct physical consequences, in 
particular logarithmically slow transport at the transition, and extensive entanglement entropy of 
the eigenstates, are consistent with recent renormalization group predictions. 

PACS numbers: 72.15.Rn, 05.30.-d, 03.75.Kk 


I. INTRODUCTION 

Experimental advances of the past decade have lead 
to the realization of isolated quantum many-body sys¬ 
tems of cold atoms and trapped ions. These systems 
have slow intrinsic time scales and are amenable to var¬ 
ious quantum-optics experimental probes, thus they are 
ideally suited for probing far-from-equilibrium quantum 
dynamics [1]. The studies of cold atoms and trapped ions 
have inspired a theoretical quest for the universal frame¬ 
work to describe non-equilibrium quantum phenomena. 
One of the central challenges that has emerged is to un¬ 
derstand the microscopic mechanisms of ergodicity and 
its breakdown in isolated quantum systems [1]. 

It has been established that there exist two distinct 
generic classes of many-body systems: the ergodic ones, 
which thermalize as a result of quantum evolution, and 
the many-body localized (MBL) [2-5], which are disor¬ 
dered and avoid thermalization via a mechanism similar 
to Anderson localization [6] in the Hilbert space. Er¬ 
godic and many-body localized systems have drastically 
different spectral and dynamical properties. Microscopic 
mechanism of thermalization in ergodic systems is called 
the Eigenstate Thermalization Hypothesis (ETH) [7, 8]. 
According to the ETH, thermalization in ergodic systems 
is manifest already in the properties of individual eigen¬ 
states, in which physical observables are thermal. In con¬ 
trast, the eigenstates of MBL systems violate the ETH 
due to the emergence of extensively many local integrals 
of motion (LIOMs) [9-13]. 

One of the invaluable tools for probing ergodicity that 
has recently emerged is the entanglement entropy. In 


ergodic systems, the eigenstates are similar to random 
vectors in the Hilbert space (modulo global conserva¬ 
tion laws, e.g. energy); therefore, the eigenstates have 
an extensive, volume-law entanglement entropy. In con¬ 
trast, the eigenstates of MBL systems can be obtained 
from product states by quasi-local unitary transforma¬ 
tions, and have low entanglement that obeys the “area 
law” [9, 14]. Further, ergodic and MBL phases are dis¬ 
tinguished by the dynamics of entanglement following a 
quantum quench from initially non-entangled states: in 
the former case, entanglement spreads ballistically and 
grows linearly in time [15, 16], while in the latter case 
the spreading is logarithmic in time [17-19]. The latter 
property has been understood from the phenomenologi¬ 
cal theory of the MBL phase based on the LIOMs [9, 10]. 

Recently, first signatures of MBL have been observed in 
systems of cold atoms [20]. In this experiment, disorder 
can be tuned in a broad range, which allows one to probe 
MBL and ergodic phase, as well as the transition be¬ 
tween them. The MBL-delocalization transition is a new 
kind of a dynamical phase transition, across which the 
eigenstates undergo a dramatic transformation as their 
entanglement entropy changes from area-law to volume- 
law. Although it is conceivable that certain many-body 
systems exhibit an intermediate non-ergodic phase, here 
we focus on the case of a direct transition between the 
MBL and ergodic phases, which is suggested by numeri¬ 
cal studies of ID systems [5, 18, 21-27], as well as by a re¬ 
cent phenomenological renormalization group study [28], 
and an effective percolation model [29] . Developing a full 
theory of this transition remains a major challenge. 

In this paper, we propose a new approach to study- 
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Figure 1. (a) The system is perturbed by a local physical op¬ 

erator V that mixes the eigenstates of the unperturbed Hamil¬ 
tonian H , leading to an effective hopping problem shown in 
the panel (b). Vertices of the lattice correspond to the un¬ 
perturbed eigenstates, and matrix elements of the operator 
V determine the hopping amplitudes. Bottom panels illus¬ 
trate the possible outcomes of a local perturbation applied to 
the Hamiltonian: all the eigenstates are mixed in the ergodic 
phase (panel c), whereas in the MBL phase the perturbation 
admixes a finite number of states (panel e). In the critical 
regime, the eigenstates are expected to be multifractal [32]. 


ing MBL, ergodicity, and the transition between them. 
We introduce a single dimensionless parameter defined 
in terms of physically measurable quantities - matrix 
elements of local operators - which provides a direct 
and sensitive probe of MBL. We formulate a criterion 
of MBL-delocalization transition in terms of this param¬ 
eter, and explore the physical properties at the critical 
point. The parameter can be viewed as a many-body 
generalization of the Thouless conductance [30] in the 
single-particle Anderson problem; the latter plays a cen¬ 
tral role in the scaling theory of localization [31]. 

The main idea of our approach, summarized in Fig. 1, 
is to study the response of the system’s eigenstates to a 
local perturbation V. We consider a lattice many-body 
system with an unperturbed Hamiltonian H and eigen¬ 
states | n), H\n) = E n \n). The problem of finding the 
spectrum and eigenstates of H + V can be viewed as a 
hopping problem on a lattice where sites correspond to 
the eigenstates | n) of unperturbed Hamiltonian iL, with 
on-site energies E' n = E n + V nn and hopping amplitudes 
V nrn between sites n,m (see Fig. 1). Thus, all the in¬ 
formation is encoded in eigenenergies E' n and in the off- 
diagonal matrix elements V nrn . Note that for local V , 
V nn < ||U||, and in general the order of energies E' n is 
different from E n . 

A key insight is that perturbing the Hamiltonian with 
a local operator V modifies the eigenstates very differ¬ 
ently in the MBL and ergodic phases [see Fig. l(a-b)]. 


In the MBL phase, the perturbation V can only affect 
the degrees of freedom within a localization length £ 
away from its support. Thus, the perturbation strongly 
admixes a given eigenstate with only a finite number of 
other eigenstates [Fig. 1(e)], as expected from the exis¬ 
tence of LIOMs in the MBL phase. In contrast, in the 
ergodic phase V hybridizes an extensive number of eigen¬ 
states, therefore its effect is highly non-local, and the 
eigenstates of H + V are completely delocalized in the 
basis of the unperturbed eigenstates. 

To detect the localization/delocalization transition, we 
ask whether hopping induced by V hybridizes the eigen¬ 
states with neighbouring values of E' n (the closest states 
in the many-body spectrum). To characterize this hy¬ 
bridization, we introduce the parameter 

Q(e,L) = In JWVU (!) 

E n +1 E n 

where we have ordered the eigenstates by their modified 
energy E ' n , and 6 = E' n /L = E n /L+0(1/L) is the energy 
density. As suggested above, the scaling of Q with system 
size L should be qualitatively different in the MBL and 
ergodic phases. 

In the MBL phase, at L £ the eigenstates close in 
energy typically have very different spatial structure due 
to their different values of LIOMs, and a local opera¬ 
tor couples them exponentially weakly in L. Therefore, 
in the MBL phase G(L) oc —L. On the other hand, in 
the ergodic phase the ETH implies that a local pertur¬ 
bation mixes the energy levels very strongly, and this 
requires G{L) oc -\-L. The above considerations suggest a 
natural criterion for the delocalization transition at the 
energy density 5 to be (G(e,L)) ~ 1. This condition im¬ 
plies that a local perturbation significantly changes the 
structure of the eigenstates, such that the LIOMs be¬ 
come non-local and the effective hopping problem enters 
the critical regime. These results will be analytically and 
numerically substantiated in Sections II and III below. It 
will be shown that (5, as a single parameter, characterizes 
the delocalization transition. 

Previous numerical studies of the delocalization tran¬ 
sition have relied on the statistics of energy levels [5], as 
well as the fluctuations of entanglement entropy [22]. In 
contrast, our approach allows us to directly probe the lo¬ 
calization properties of the eigenstates through the ma¬ 
trix elements of local operators. Unlike entanglement 
entropy and level statistics, these are directly related 
to correlation functions, and therefore are, in principle, 
measurable. Moreover, our approach yields further con¬ 
sequences for the properties of the system at the critical 
point. Those include the extensive entanglement in the 
eigenstates, and logarithmic in time transport of con¬ 
served quantities and entanglement propagation. These 
properties are qualitatively consistent with the predic¬ 
tions of recent phenomenological renormalization group 
studies [28, 29]. 

The rest of the paper is organized as follows. In Sec¬ 
tion II we provide analytical derivation for the scaling 
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of Q with energy density and system size in the MBL 
and ergodic phases, and formulate the delocalization cri¬ 
terion. Using exact diagonalization, in Section III we 
then test our approach numerically on the random-field 
XXZ spin chain, where the transition between MBL and 
ergodic phases can be driven by changing the disorder 
strength. We demonstrate that the G(L) ~ 1 crite¬ 
rion gives a location of the transition that is in good 
agreement with the estimates obtained using other meth¬ 
ods [24]. Furthermore, we show that our approach has 
a good energy-density resolution, and apply it to map 
out the phase diagram as a function of disorder and en¬ 
ergy density, demonstrating the presence of the many- 
body mobility edge in the random-field XXZ model. In 
Section IV we address the physical consequences for the 
MBL transition that follow from our analysis. In par¬ 
ticular, we find the logarithmic in time propagation of 
entanglement and particle number fluctuations, in agree¬ 
ment with recent works [28, 29]. We numerically verify 
those predictions using exact diagonalization and time- 
evolving block-decimation simulations of a global quench 
in the random-field XXZ model. Our conclusions and a 
discussion of future directions is presented in Section V. 


II. STATISTICS OF MATRIX ELEMENTS AND 
THE SCALING OF Q: ANALYTIC 
CONSIDERATIONS 

In this Section we discuss the distribution of the pa¬ 
rameter Q introduced in Eq. (1) and its scaling in the 
MBL and ergodic phases. In the former case, the be¬ 
haviour of G can be understood by invoking the ETH. 
In the latter, our considerations are based on the picture 
of LIOMs in MBL phase. In the subsequent Section, we 
present data from exact diagonalization which are fully 
consistent with these analytical arguments. 


A. Ergodic phase: matrix elements from the ETH 

In the ergodic phase, the eigenstates are expected to 
be similar to random vectors in the Hilbert space (this 
analogy is qualitative due to locality of interactions and 
energy conservation). 

Consider two eigenstates |n), | m) near the middle of 
the many-body band. They satisfy the orthogonality 
constraint (m\n) = 0, n V m. Intuitively, acting on 
the state | n) with a local perturbation, V\ n), gives an¬ 
other random vector in the Hilbert space, which is no 
longer orthogonal to | m). Thus, it is expected that the 
off-diagonal matrix element V nrn = (m\V\n) is a scalar 
product of two random vectors, | m) and V\n). Denot¬ 
ing the Hilbert space dimension by D, the matrix ele¬ 
ment can be estimated as V nrn ~ 1/VU. For a quantum 
many-body system, in general V ~ exp(sL), where L 
is the number of lattice sites or spins, and s is statis¬ 
tical entropy density. Since the level spacing scales as 


the inverse Hilbert space dimension, 5 ~ 1 /D, we have 
S < Vim, and a local perturbation strongly mixes the 
nearby eigenstates. Our parameter Q , which in the lead¬ 
ing order is given by G(L) ~ In (V nm /S) ~ 1/2 In V cx L, 
thus grows linearly with system size. 

The intuitive argument above can be extended to 
states with an arbitrary energy density using Sred- 
nicki’s [33] ansatz for the off-diagonal matrix elements 
of local operators: 

v nm = e~ s ^ E ^' 2 f(E n , E m )R nm , (2) 


where S(E) is the statistical entropy at energy E = {E n + 
Em)/2 , R n m is a random number of order one, and / is 
a smooth function which can be linked to the thermally- 
averaged response of the system to the perturbation V. 
Since S oc e~ s ^ E \ we obtain that 


G(L) = In 


Vr 


n,m S(E,L) 

~— oc - oc 


s(e) 


L, 


(3) 


for states with an arbitrary energy density e = E/L , 
where s(e) is the corresponding entropy density. It is 
worth noting that the off-diagonal matrix elements in the 
ergodic phase have been studied numerically in Ref. [34] , 
where the scaling V nm oc l/VV was confirmed. 


B. MBL phase: matrix elements from local 
integrals of motion 

To understand the properties of the off-diagonal ma¬ 
trix elements V nm and Q deep in the MBL phase we rely 
on the picture of LIOMs [9-12]. It was argued that the 
key property of a system that shows many-body local¬ 
ization at infinite temperature is that its eigenstates can 
be connected to product states by a quasi-local unitary 
transformation [/, which entangles the remote degrees of 
freedom only exponentially weakly. 

For concreteness, let us consider an interacting spin- 
1/2 model [e.g., the XXZ model introduced in Eq. (8) 
below]. In this case, it is convenient to choose the nat¬ 
ural “up-down” basis, i.e., the basis vectors are eigen¬ 
states of all erf operators. Then, in the MBL phase 
there exists a quasi-local unitary operator U (quasi-local 
means that this operator only creates short-range entan¬ 
glement), which diagonalizes the Hamiltonian in this ba¬ 
sis, WHU = Hdiag- The unitary U defines a mapping 
from the physical spin operators, erf, to “effective spins” 
rf via 


r? = UafUl (4) 

The operators rf form a complete set of mutually com¬ 
muting, quasi-local integrals of motion in the MBL phase. 
In the r-basis the eigenstates therefore simply become 
“up-down” states of the effective spins. 

In order to analyze the effect of a local perturbation, 
we expand the corresponding local operator V (which for 
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simplicity is chosen to act on the first few spins in the 
chain) over the complete basis of operators 1 t^ 2 ...t/* I/ , 
where E {0, d=, Here r 2 ° = 1 denotes the iden¬ 
tity operator on the site i, and =b denote the usual rais¬ 
ing/lowering operators. We represent V as a sum over 
operators V r of range r, which affect only the effective 
spins 1,..., r: 

v = v r - v; ai - a -rp...r“-, (5) 

r ai,...,a r 

where the last a r ^ 0. The off-diagonal matrix element 
of V between the eigenstates which differ by ~ L flips of 
effective spins is equal to where ~ L coefficients 

are ai = ± such that the operator 1 ... r^ r maps one 
configuration into another. For example, if |ni) 
and 1 77 , 2 ) the matrix element (n2\V\ni) = V 1 *" 

gives us a coefficient in front of term in the 

expansion (5). Thus the matrix elements of V probe the 
nature of the spatial decay of the coefficients V^ 1 '" 01, , 
which, in turn, encode the properties of the unitary V 
(in particular, they control the amount of entanglement 
generated by U). 

Deep in the MBL phase, the LIOMs are robust un¬ 
der local perturbations. Typically, V hybridizes a given 
eigenstate only with states which have the same inte¬ 
grals of motion distance x < £ (£ being the localization 
length) away from where V acts. Two generic eigenstates 
typically have different integrals of motion everywhere. 
Hence, in order for a local perturbation V to couple these 
states, it has to flip the state of effective spins across 
the entire volume of the system. This requires tunneling 
an “excitation” through the entire system, which would 
leave behind a trace of spin flips [35]. In the localized 
phase, such tunneling is exponentially suppressed; as a 
consequence, the probability of V to couple such states 
decays exponentially with system size. Thus, typically 
the ratio V nrn /d behaves as c ~ kL , and the parameter Q 
decays linearly with the system size: 

G(L) oc —kL. (6) 

Below we study Q (L) numerically for the XXZ spin chain, 
finding the behaviour consistent with the above formula. 
We also find that G(L) has a broad, nearly normal distri¬ 
bution, similar to the distribution of the wave function 
tails in the single-particle Anderson insulator [36]. 

Note, that Eq. (6) indicates that matrix elements of 
a local perturbation are similar to a random ultramet¬ 
ric matrix [37]. Indeed, the scaling of Q in Eq. (6) 
implies that a local perturbation which changes values 
of the LIOM up to a distance x away from its cen¬ 
ter [there are 2 X ~ 1 such eigenstates], has a matrix ele¬ 
ment V(x) oc e -(^+ s )^ where s is the entropy density. 
For k > 0 the eigenstates are localized, and the transi¬ 
tion occurs at k = 0 [37], which is similar to our crite¬ 
rion. The problem of Anderson localization for random 
ultrametric matrices has been studied analytically using 
strong-disorder renormalization group [38], and, in par¬ 
ticular, multifractal exponents have been computed [37]. 


Note, however, the im- portant difference of the ensem¬ 
ble of ultrametric matrices considered in the literature, is 
that the matrix elements there obey the gaussian, rather 
than the log-normal distribution. We reserve the explo¬ 
ration of this analogy to further studies. 

We note that the behavior of matrix elements in the 
MBL phase was used [39] to study the properties of a 
local spectral function. In addition, the spatial decay of 
the diagonal part of local operators given by terms with 
all a.i E {0 ,z} in Eq. (5) was studied in Ref. [13], and 
was also found to decay exponentially in the MBL phase. 

C. A criterion for the transition 

Using the predicted linear scaling of the parameter Q 
with system size in the ergodic and MBL phases, it is 
natural to propose a phenomenological criterion for the 
delocalization transition as the absence of scaling of Q 
with L : 

G{L) = Sc (7) 

Below, we will use this criterion to determine the loca¬ 
tion of MBL-ergodic transition as a function of disorder 
strength and energy density in exact diagonalization. 

Physically, the above criterion implies that LIOMs, 
characteristic of the MBL phase, become unstable and 
are expected to delocalize at the transition. Indeed, any 
local perturbation strongly couples a given eigenstate 
with an adjacent eigenstate (in the reshuffled spectrum 
E[ J, which has a different spatial structure, and there¬ 
fore different LIOMs. In this case, the new eigenstates 
are superpositions of states with different values of the 
LIOMs across the entire system. Thus, a local perturba¬ 
tion modifies the integrals of motion even far away from 
where it is applied, which implies that integrals of motion 
necessarily become non-local. 

We have argued that the above condition implies the 
non-locality of integrals of motion. It is natural to ask 
whether the criterion (7) indeed implies delocalization in 
the effective hopping problem we introduced in Section I 
above, and what are the properties of the system at this 
delocalization transition. That is, whether the participa¬ 
tion ratio of new eigenstates in the basis of unperturbed 
eigenstates tends to infinity as L —>> 00 . One cannot an¬ 
swer this without additional information regarding the 
distribution of matrix elements. We have studied several 
quantities which probe the delocalization in the hopping 
problem directly, and indeed found that the above crite¬ 
rion is consistent with the delocalization in the hopping 
problem. 

III. NUMERICAL RESULTS 

In this Section we test the predictions for the distribu¬ 
tion and scaling of Q using numerical (exact) diagonaliza¬ 
tion of the one-dimensional random-field XXZ spin-1/2 
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Figure 2. Distribution of Q = log (\V nn +i\/5) across the MBL transition displays qualitatively different scaling with system 
size, (a) For weak disorder (W = 0.5),when the system is in the ergodic phase, Q increases with system size, and the distribution 
shifts to the right, (b) At the MBL transition (W = 3.6), the distribution broadens but does not move, (c) In the MBL phase 
(W = 5), G becomes smaller for larger systems, and the shape of the distribution is approximately gaussian. 


chain. The Hamiltonian of this model is given by: 

L L L 

H=J X 53 K + a~ +1 +a~ <jf +1 )+-J z 53 ofof +1 +53 V, 2 , 

i= 1 i= 1 i 

(8) 

where we use periodic boundary conditions and identify 
<jl+i = a i- Below we set J x = 1. The random z- field hi 
on each site is drawn from a box distribution [—FF; W], 
and a ± = (a x ± icr v )/2. At J z = 0, this model can be 
mapped onto free spinless fermions moving in a disorder 
potential via the Jordan-Wigner transformation; in this 
case, all eigenstates are localized for arbitrary disorder 
strength W. Introducing J z ^ 0 is equivalent to turning 
on the nearest-neighbour interactions between fermions. 

Previous work [5] has demonstrated that the model (8) 
exhibits the MBL phase at strong disorder (W > W c ), 
and the ergodic phase at weak disorder (W < W c ). The 
delocalization transition takes place at some critical dis¬ 
order strength W c . Several numerical studies [5, 21, 23- 
27, 40, 41] have estimated the location of the transi¬ 
tion in this model using various probes, including the 
statistics of eigenenergies, entanglement entropy and its 
fluctuations, participation ratios and different dynami¬ 
cal probes. In particular, state-of-the-art exact diagonal- 
ization performed on systems up to L = 22 spins [24] 
has robustly determined the location of the transition at 
W c « 3.6. 

Here we use exact diagonalization to numerically study 
the statistics of the off-diagonal matrix elements and pa¬ 
rameter G for the XXZ spin chain (8) across the MBL- 
delocalization transition. In order to establish the uni¬ 
versality of this approach, we considered the response of 
the spin chain to two different local perturbations: 

Vi =af, V 2 = aj aj +1 + aj af +1 . (9) 

The interaction strength was fixed to be = 1, simi¬ 
lar to previous studies. All the data were obtained for 
chains with periodic boundary conditions constraining 
the Hilbert space to states with zero total spin, such that 


the dimension is 'D(L) = ( L / 2 )* Disorder averaging was 

performed over at least 10 4 , 5000, 4000, 2000 and 200 
realizations for L = 8,10,12,14,16 spins. 

A. Distribution of G(L) in the middle of the band 

First we study the disorder-averaged distribution of G 
over V(L)/L eigenstates in the middle of the band. As 
discussed in the introduction, we perturb the eigenspec- 
trum by the diagonal matrix elements of the chosen op¬ 
erator, E' n = E n + (n|Vi ? 2 1^)• Perturbed eigenenergies 
{E' n } are sorted, and we collect the statistics for G(L) as 
defined in Eq. (1). 

The behaviour of the distribution function of G(L) is 
illustrated in Fig. 2. In this case, the perturbation opera¬ 
tor was chosen to be V\ (we found a qualitatively similar 
behaviour for V 2 ). We observe that the distribution of 
G exhibits a qualitative change across the MBL transi¬ 
tion. In the ergodic phase [Fig. 2(a)], the average value 
of ( G{L )) grows with the system size, signalling that a lo¬ 
cal perturbation is more and more likely to hybridize the 
nearby states - a signature of delocalization. If one scales 
the distribution by the square root of the Hilbert size di¬ 
mension, (Q(L))/\/V, the distributions for different sys¬ 
tem sizes indeed approximately collapse, confirming that 
G scales as predicted by the ETH. 

In contrast, in the MBL phase [Fig. 2(c)], the averaged 
value of ( G) decreases as system size is increased; this 
reflects the fact that a local perturbation is less and less 
likely to hybridize the nearby states, and therefore the 
eigenstates and the LIOMs remain robust. Moreover, 
the distribution of G becomes broader at larger L, and 
is approximately normal. We found that the distribution 
becomes very close to normal at large system sizes and 
strong disorder. 

Finally, at intermediate disorder W = 3.6, the aver¬ 
age (G(L)) remains nearly independent of system size. 
Hence it is natural to identify this point with the MBL- 
delocalization transition. Below we will use the energy- 








6 


resolved average value of Q to identify the transition point 
and the mobility edge. Note that although the average 
value of (G(L)) remains constant at the transition, the 
distribution broadens as the system size is increased. De¬ 
tailed investigations revealed that this effect originates 
largely from sample-to-sample fluctuations, rather than 
from the state-to-state fluctuations in the same disorder 
realization. 

B. Average energy-resolved (Q(e,L)) 

In order to obtain the location of the MBL transition 
more precisely, we now study the average value (G(s, L )} 
as a function of the dimensionless energy density and 
system size. The dimensionless energy density 5 (referred 
to as “energy density” in the following, for simplicity) is 
defined as £ — (E -E^min)/ (E max E m i n ), where E/ m j n / max 
are the energies of ground state (highest excited state) of 
our system. 

Figure 3 shows the system-size dependence of 
(G(£ c , L)) for fixed s c = 0.45, which is the energy density 
at which the delocalized phase is most robust. Similar to 
the behaviour already observed for the distribution of £?, 
we see that the behaviour of averaged (G(£ C i L)) is quali¬ 
tatively different at weak and strong disorder. At W < 3 
we have d(G(s c , L))/dL > 0, and the second derivative 
appears to be positive, signalling that larger systems be¬ 
come more and more thermal. At strong disorder W > 4, 
(G(£a L)) behaves according to Eq.(6), as expected in the 
MBL phase. From Fig. 3 we identify the critical value of 
disorder 

W c = 3.6 ±0.15 

using our criterion for the MBL transition, Eq. (7). This 


agrees with the previous findings of Refs. [24]. 

Using the same delocalization criteria, we can map out 
the many-body mobility edge in the random-field XXZ 
model. We define the many-body mobility edge to be at 
the energy density e(W) where (G(s c ,L)) is independent 
of system size. Results of this calculation are shown in 
Fig. 4. We also notice the intrinsic asymmetry of the 
mobility edge: it is shifted towards the ground state and 
centered slightly below £ = 1/2 [24, 40, 42]. This asym¬ 
metry becomes stronger with increasing the interaction 
strength (not shown), which agrees with the recent argu¬ 
ments of Ref. [42] . 


IV. ENTANGLEMENT AND DYNAMICS AT 
THE DELOCALIZATION TRANSITION 

Next, we explore the physical properties of the system 
at the transition, defined using the above criterion. In 
particular, we will be interested in understanding the en¬ 
tanglement properties of the eigenstates, as well as the 
dynamical properties - spreading of entanglement and 
transport of conserved quantities. We discuss qualitative 
expectations and support them with numerical studies. 

In order to qualitatively understand the entanglement 
properties of the eigenstates, we consider a system of size 
2 L and initially disconnect it across the middle link into 
left and right subsystems. The eigenstates of such a sys¬ 
tem are simply product states of the eigenstates in the 
left and right parts, \ti)l 0 \^)r, n,m = 1 ,..., V. Restor¬ 
ing the coupling between L and R systems corresponds 
to a local perturbation Vlr = ^2 a Wl 0 V a R , where the 
sum includes a small finite number of terms (three for 
the case of the XXZ spin chain), and F q l, VaR act only 
on the degrees of freedom in the left/right part. The 



Figure 3. Scaling of C/(e c , L) in the middle of the band with 
the system size for different disorders W = 1.5,...,5. Value 
of disorder is shown on the right of each curve. From here the 
critical disorder strength is determined as W c = 3.6 ± 0.15. 



2.5 3 3.5 4 
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Figure 4. Many-body mobility edge e(W) as a function 
of disorder. Blue (red) color indicates regions where (g(s : L)) 
decreases (grows) with L. Yellow regions correspond to points 
where we cannot determine the behavior due to errorbars. 
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problem of finding the eigenstates reduces to a hopping 
problem, where the sites are product states \u)l 0 I'm) r, 
and hopping amplitudes are set by the operator Vlr- 
This problem bears many similarities with the case of a 
local perturbation acting on a single system, discussed 
above. 

Previously in Section II we established that the MBL 
transition corresponds to the critical regime of the ef¬ 
fective hopping problem. In this case, we expect that 
the eigenstates of the connected system 1 1) are mul¬ 
tifractal when expressed in the product basis 1 1) = 
E n ,m^nml n )^ ® | m) r, with disorder-averaged partici¬ 
pation ratios satisfying a relation P q = m |A^ m | 2gr oc 
V~ 2Tq . This is consistent with recent work [43] which 
studied dynamical signatures of the fractality near the 
MBL transition. We have tested and confirmed the mul- 
tifractality of the hopping problem numerically (details 
will be presented elsewhere). Although we are not aware 
of a direct relation between entanglement entropy and 
participation ratios in the many-body case (for the single¬ 
particle problem such a relation does exist, see Ref. [32]), 
it is natural to expect that the entanglement entropy of 
half the system will be extensive, but sub-thermal: 

*^ent (-^ 7 ) ^ asL, OL <C 1, 

where sL is the maximal (thermal entropy) of the sub¬ 
system of size L. This agrees with the predictions of 
Ref. [28]. We note that the above relation is expected 
to hold for entanglement between two regions of approx¬ 
imately equal size. If, on the other hand, one of the 
subsystems is much smaller than the other one, its en¬ 
tanglement entropy should approach the thermal value, 
in accordance with the general arguments based on the 
strong subadditivity of entanglement [44]. 

Next, we discuss the dynamics at the transition. As¬ 
suming the analogy with the problem of ultrametric ran¬ 
dom matrices, described above (although it should be 
kept in mind that this analogy is almost certainly incom¬ 
plete - in particular, in our case the off-diagonal matrix 
elements follow a log-normal distribution rather than a 
random distribution), we expect that at the transition 
the typical relaxation time scale for a system of size L is 
exponentially long in L, 

^~m =eSL ’ (10) 

and therefore the transport at the transition will be 
logarithmically slow, consistent with the predictions of 
Refs. [28, 29]. Note that this is a time scale for both 
dephasing , and the on-shell dynamics (transport) at the 
MBL transition. For example, let us consider a setup 
in which the left and right subsystems are initially dis¬ 
connected, i.e., are prepared in an eigenstate |^(0)) = 
\uo)l \mo)R, and the coupling between them is turned 
on instantaneously at t = 0. In the MBL phase after 
a quench, due to dephasing, the entanglement entropy 
increases logarithmically, and this growth is unbounded 


in an infinite system; but there is no transport (except 
for rare resonances). At the transition, however, there is 
transport of the z—component of spin polarization, which 
can be characterized by its fluctuations in the left subsys¬ 
tem as a function of time. Hence, both fluctuations of the 
magnetization and entanglement growth are unbounded 
in the thermodynamic limit and increase logarithmically 
in time, 

Sent(t) ~ In t, (S 2 z p - (S Z , L ) 2 ~ In t, (11) 

where S z is the total z—projection of spin in the left 
part. We expect that these laws will be generic and will 
hold for a class of initial states, including the case when 
the system is prepared in a product state (rather than in 
a product state of eigenstates of the left and right parts). 

We have studied the entanglement growth and spin 
transport numerically to verify Eq. (11). We considered a 
global quench in which the system is initially prepared in 
a product state (random up-down configuration of spins), 
and time evolved with the Hamiltonian in Eq. (8). Us¬ 
ing exact diagonalization, we have been able to access 
long-time dynamics in systems of up to L = 18 spins. 



Figure 5. Dynamics of entanglement entropy (a) and parti¬ 
cle number fluctuations (b) across the MBL transition. The 
system was initialized in a product state where each spins is 
up or down. For entanglement entropy, the transition cor¬ 
responds to faster-than logarithmic growth of entanglement. 
The particle number fluctuations grow logarithmically in the 
vicinity of the transition. 











In addition, we have performed time-evolving block dec¬ 
imation [45] (TEBD) simulations which allow access to 
much larger spin chains, but are restricted to moderate 
times t ~ 10 1 2 [46]. For the TEBD algorithm we used 
the Trotter decomposition with time step At = 0.01 
and a maximum bond dimension x = 2000. The re¬ 
sults on the entanglement growth and the fluctuations of 
the spin polarization, obtained using TEBD for L = 24 
spins, are illustrated in Fig. 5. Both quantities spread 
logarithmically in time. We note that slow transport 
opens an additional channel for entanglement growth, in 
addition to dephasing [18], and as a result, the entan¬ 
glement S ent (t) near the transition does not display a 
characteristic “knee” which corresponds to the crossover 
between transport-dominated and dephasing-dominated 
growth deeper in the MBL phase. 

V. SUMMARY 

We presented a new approach for studying the MBL 
and delocalization transition based on the response of 
system’s eigenstates to a local perturbation. We char¬ 
acterized the effect of a local perturbation by a dimen¬ 
sionless parameter Q , defined via the ratio of a matrix 
element of local operator to the level spacing [Eq. (1)]. 
This parameter is reminiscent of the Thouless conduc¬ 
tance in the single-particle localization problem, which 
describes the response of energy levels to changing the 
boundary conditions. Similar to the single-particle case, 
we argued that the transition can be identified by the 
relation Q{L) ~ G c ~ 1. However, it should be noted 
that, unlike the Thouless conductance, the parameter Q 
is not directly related to the physical conductance of a 
many-body system. 

We have verified our approach by studying the MBL 
and delocalization transition in the random-field XXZ 
model. We obtained an estimate for the critical disorder 
strength at the transition which agrees very well with 
previous numerical studies based on entanglement and 
energy level statistics. Further, we mapped out the mo¬ 
bility edge in this model, demonstrating that the method 
has a good energy resolution, while at the same time re¬ 
quiring little beyond exact diagonalization of the Hamil¬ 
tonian. 

In addition to being an efficient tool to detect the MBL 
transition, our criterion gives insights into its microscopic 
details. In particular, our criterion implies that at the 


MBL transition, both transport of conserved quantities 
(such as spin polarization), as well as spreading of en¬ 
tanglement, are logarithmically slow. Such dynamics at 
the MBL transition was predicted by recent phenomeno¬ 
logical RG studies of the MBL transition [28, 29]. Using 
exact diagonalization and TEBD calculations, we have 
computed the spreading of entanglement and transport 
in large systems, confirming the logarithmic in time be¬ 
haviour for both quantities. 

Finally, we mention some further directions opened 
by this study. First, we expect that this method will 
be useful for studying MBL and ergodicity breaking 
in other contexts, for example in the translationally- 
invariant models which were conjectured to break ergod¬ 
icity [47-53]. Second, our results provide a natural start¬ 
ing point for developing a microscopic renormalization 
group (RG) procedure for the MBL transition. In spirit, 
the microscopic parameter Q is similar to the phenomeno¬ 
logical parameter (ratio of the “entangling rate” to the 
level spacing), introduced in a recent RG study [28], while 
having the advantage of being numerically measurable. 
Hence, one could potentially use the “physical” distri¬ 
bution of Q, obtained from exact diagonalization, as an 
input parameter, and make use of the RG procedure to 
reach larger length scales. This would complement the 
existing phenomenological RG studies and is left for fu¬ 
ture work. Another interesting direction would be to 
further explore the consequences of the LIOMs in the 
MBL phase. As we discussed above, our results hint at 
an underlying ultrametric structure, which could allow 
one to introduce a random matrix ensemble from which 
the properties of the MBL phase and the delocalization 
transition could be analytically computed. 
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